close all; 
clear clc; 
Ip0=imread('E:\tuxiang_shiyan\zhiliao\Sample_Pictures\Desert2.jpg');
Ip=rgb2gray(Ip0);
[ma,na]  = size(Ip);
%最大类方间差
a = Ip;                                 
[m,n]=size(a);                           % 图像大小
max_g = 0;
%t = 180;
for t = 0 : 255
    N0 = 0;
    N1 = 0;
    u0 = 0.0;
    u1 = 0.0;
    a = double(a);                       % 强制转化,否则平均灰度值将无
                                         % 法计算
    for i = 1:m   
        for j = 1:n   
            if a(i,j) <= t   
                N0 = N0 + 1;         % 小于阈值的像素点数目
                u0 = u0 + a(i,j);    % 目标总灰度
            else   
                N1 = N1 + 1;         % 大于阈值的像素点数目
                u1 = u1 + a(i,j);    % 背景总灰度
            end   
        end   
    end
    w0 = N0 / (m*n);                    % 初始化
    w1 = N1 / (m*n);
    u0 = u0 / N0;
    u1 = u1 / N1;
    u = w0 * u0 + w1 * u1;
    g = w0 * (u0 - u)^2 + w1 * (u1 - u)^2;
    if (g >= max_g)
        max_g = g;
        T = t;
    end
end
% 二值化
for i=1:m   
    for j=1:n   
        if a(i,j) < T   
            a(i,j) = 0;   
        else   
            a(i,j) = 255;   
        end   
    end   
end
figure(2);
subplot(1,2,1);
imshow(mat2gray(Ip0));
title('原图');
subplot(1,2,2);
imshow(mat2gray(a));
title('分割结果');